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We study numerically the time evolution of the transport properties of layered superconductors after 
different preparations. We show that, in accordance with recent experiments in Bi2Sr2CaCu20g 
performed in the second peak region of the phase diagram (Portier el at, 2001), the relaxation 
strongly depends on the initial conditions and is extremely slow. We investigate the dependence on 
the pinning center density and the perturbation applied. We compare the measurements to recent 
findings in tapped granular matter and we interpret our results with a rather simple picture. 



The behavior of the flux-line array in type II super- 
conductors is determined by two competing interactions, 
the vortex-vortex repulsion that favors order in the form 
of a hexagonal lattice, and the attraction between vor- 
tices and randomly placed pinning centers due to mate- 
rial defects that, as thermal fluctuations, favors disorder. 
This competition generates a very rich static phase di- 
agram [|]. More recently, the evolution of the spatial 
ordering of vortices driven out of equilibrium by, e.g., an 
external current has been analyzed analytically 1^, nu- 
merically and experimentally jj]. 

However, even in the absence of an external force, vor- 
tex systems show very rich nonequilibrium phenomena. 
Indeed, a hallmark of vortex dynamics is its history de- 
pendence. A vortex structure prepared by zero-field cool- 
ing (ZFC), in which the magnetic field is applied right 
after crossing the transition temperature, is expected to 
be more ordered than the one prepared by field cooling 
(FC), in which the magnetic field is applied before going 
through the transition |5|-|^] (see, however, [||). Subse- 
quently, thermal fluctuations allow for the rearrangement 
of the initial configurations. Studies of the evolution af- 
ter FC and ZFC preparations of low- Tc and high-Tc type 
II superconductors have been performed with a variety of 
techniques that we non exhaustively summarize below. 

The initial current ramp in the IV characteristic of a 
FC 2H-NbSe2 sample has a large hysteretic critical value, 
suggesting that vortices are strongly pinned ^,0]. Once 
the fiux lines are depinned, and one applies a subse- 
quent ramp, the critical current takes a smaller value. 
In ZFC samples the critical current in both ramps is 
small The presence of metastability and hysteresis 
depends on the speed of the current ramp imposed to 
the system The ac response and complex resistiv- 
ity in 2H-NbSe2 also show that the FC state is strongly 
pinned and disordered while the ZFC state is not 0. 
Similar effects have been exhibited in high Tc supercon- 
dutors. For instance, Josephson plasma resonance ex- 
periments in Bi2Sr2CaCu208 (BSCCO) suggest that the 



resonant field in FC samples is stationary while in ZFC 
samples it approaches the FC value asymptotically [ p^ . 
Still, thermomagnetic history effects in the solid vortex 
phase of pure twinned YBa2Cu307 single crystals have 
been recently observed in ac susceptibility measurements 
|]ll| . On the numerical side, Olson et al studied how the 
fixed drive induced velocity (voltage response) depends 
upon time after different sample preparations, close to 
the order (3D) - disorder (2D) transition focusing 
on the effect of superheating or supercooling the ordered 
and disordered phases |p^ . 

In this Letter we concentrate on a recent study of 
the long-time transport properties of BSCCO monocrys- 
tals after different preparations |lj]. The experimental 
protocol is as follows. The same working conditions, 
given by a temperature T = A.bK and a magnetic field, 
fioH — 1.5T, perpendicular to the c-axis are initially 
attained via a FC or a ZFC procedure. These values 
fall in the second peak region After a waiting-time 
tiu — 30min a triangular pulse of current of duration 
Tp ~ 10/is — 500/is is applied and an IV characteristic 
is recorded. During a second waiting period of the same 
duration 1^ the sample is let freely evolve driven only 
by thermal fiuctuations. The same external current is 
then applied and a second IV characteristic is recorded. 
This procedure is repeated so on and so forth. The first 
observation is that the shape of the IV loop changes as 
time elapses. In order to quantify its change, Portier 
et al chose a threshold tension, ^th, defined the thresh- 
old current, Jth, as the corresponding intensity in the IV 
characteristics, and monitored the dependence of Ijh on 
the discrete times t„ — n{tw -f- Tp), n = 0,1,.... /th 
slowly relaxes in time, in a way that depends on the his- 
tory of the sample. In the FC sample, after a seemingly 
static period that lasts until t w 10* s, 7™ relaxes ~ 25% 
of its initial value over one decade in time, with a loga- 
rithmic decay. In the ZFC sample instead, /th smoothly 
increases by less than a tenth of its initial value over 6 
decades. The asymptotic value of /th is the same for 
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both coohng procedures. The samples evolves in very 
long time scales, of the order of days. 

The microscopic interpretation of experiments that 
probe the flux motion via the measurement of transport 
properties is not always straightforward. We reproduced 
the latter experiment with numerical simulations using 
a model of pancake vortices in a layered superconductor 
that successfully captures many of the observed effects in 
type II superconductors [ p^lp^ -p^. For adequately cho- 
sen sets of parameters our results are in qualitative ac- 
cord with the measurements of Portier et al . The ad- 
vantage of using numerical simulations is twofold. First, 
it allows us to simply explore the effect of the microscopic 
parameters in the model. Second, it allows us to grasp 
what vortices are actually doing by direct visualization. 

The model takes into account the long-range magnetic 
interactions between all pancakes. Since a flux line is 
essentially massless [Q we use overdamped Langevin dy- 
namics. The equation of motion for a pancake located at 
a position = (r^, z^) = (xi,yi, Zi), (z = c), is 



(1) 



where — \ri — Vj] and Zij = \zi — Zj\ are the in- 
plane and inter-plane distances between pancakes i and 
J. ?7 is the friction coefficient. The driving force per unit 
length due to an in-plane current J acting on a vortex 
is F = $0 J X z, with $o the quantum of magnetic flux. 
Quenched random point disorder is modeled by pinning 
centers that occupy uncorrelated random positions, taken 
from a uniform distribution, on each layer. They exert 
an attractive force on the vortices, 
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where Ap measures the amplitude of this force, Op is the 
pinning range and r^p — jr^ — rp| is the in-plane distance 
between vortex i and a pinning site at Rp = (rp, Zj). The 
magnetic interaction between pancakes is 
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with R = \/ z"^ + r^, A — 2A^/(i the 2D thin-film screen- 
ing length, d the inter-layer spacing and A|| the in- 
plane penetration length. One can vary the relative 
strength of the intcrlayer coupling by tuning the prefactor 
s [|l|,|l|,|l8|. The model in Eqs. (|) and (|) is vahd in the 
limit d ^ A| I <C A |^ . We normalize length scales by A| | , 
energy scales by Ay = $Q/47r^A and time by r = rjX^^/Ay, 
that for BSCCO yields r ~ 1/is. We consider iV„ vortices 
and Np pinning centers per layer in Ni rectangular layers 
of size Lx X Ly {Lx = \/3/2Ly). The normalized vortex 



density is Uy = B/^q = Ny/{LxLy). Following ||T^ we 
consider Uy = 0.12, Ap/Ay = 0.2, ap = 0.1, d = 0.01 
and A = 200 which correspond to BSCCO |§. We use 
Ny = 25 vortices in up to Ni = 4 layers. The value of Qp 
as compared to Hp is of first importance and will be dis- 
cussed later on. We use periodic boundary condition in 
all directions and the periodic long-range in-plane inter- 
action is evaluated with an exact and fast converging sum 
]l7t . All the simulations are at T = 0, with a time-step 
0.01, and averaged over 90 realizations. 

We mimic the external current applied with an in- 
creasing (in steps) force I parallel to the x axis. It re- 
sults in an average motion of the pancakes, with velocity 
y = [{Vxim = {NyNl)-^Y.^[{dx^/dt)]. (The angular 
and square brackets represent an average over different 
realizations of the dynamics and disorder, respectively.) 
We match the notations in Ref. jlj] by calling /th the 
value of the applied force which leads to an average speed 
equal to a previously chosen threshold 1th- We have 
checked that the results do not depend qualitatively on 
the precise value of Vth, if e.g. Kh < 0.2 for /max — 1-4 
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FIG. 1. Current-voltage curves at times = (open cir- 
cles) and t2 = 20 (filled circles) for a monolayer of pancake 
vortices with rip = 0.12, = Tp = 14 and /max = 1-2. 

Experimentally, vortices are created on the surface of 
the sample when the magnetic field is turned on. It hap- 
pens during cooling for FC samples and later for ZFC 
samples. Afterwards, they penetrate the sample and new 
ones are created on the border. This procedure is diffi- 
cult to implement in this model. We differentiate the FC 
and ZFC initial states by choosing different initial repar- 
titions of vortices. For the FC samples, we choose an 
initial condition with each vortex pinned to a defect, in 
which the distance to the nearest pinning center is less 
than the pinning range ap. To reproduce the ZFC sam- 
ples, we choose an ordered hexagonal lattice This 
procedure neglects the critical state and it is expected to 
describe the behavior in the bulk. 

In Fig. we display the first two IV loops for a 2D sys- 
tem that starts from a pinned (FC) or hexagonal (ZFC) 
initial condition. As expected [pyi9|, the first ramp shows 
hysteresis and a large critical current for the FC sample. 
In the second ramp the loop moves and it is partially 
closed. The ZFC sample has a much lower critical cur- 
rent (that may even vanish). The pulse we used is as 
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fast as the fastest used experimentally, Tp = 14 in real 
time; consequently, the critical currents are very weak. 
Moreover, since the waiting-time between two successive 
pulses is here t^, Tp, while in t^, ^ Tp, we are 
perturbing the system in a much stronger manner than 
done experimentally. By the end of this Letter we shall 
study the dynamics under a less invasive probe. 

The time-dependence of the IV loop is quantified in 
Fig. ||-lcft where we show the evolution of /th against 
tn for two Up. The curves correspond to pinned (open 
circles) and hexagonal (filled squares) initial conditions. 
The first two points in the lower curves for the pinned and 
hexagonal initial conditions correspond to the IV curves 
in Fig. There is an initial fast motion that translates 
into a rather fast drop in I^^ and a fast increase of /th*^. 
Later, the curves approach their common asymptote in 
a slower than exponential manner. The qualitative be- 
havior of these curves resembles the ones in ||lj] though 
with much shorter relaxation times. 
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FIG. 2. Left: The threshold current /^(tn) against time 
in a monolayer with rip = 0.12, 0.25 (upper to lower asymp- 
totic value). Open symbols: pinned, FC; filled symbols: lat- 
tice, ZFC. /max = 1.4. Right: /^(tn) for /max = 1.4, Tp = 14, 



= 6 (above) and /max = 1.6, Tp = 26, 
(open symbols) and ZFC (filled symbols). 



= 4 (below). FC 
i-o = 0-25. 



A closer inspection of the model and procedure sug- 
gests that the precise shape of these curves may depend 
on several parameters like Up, the form of the ramp and 
the choice of the sequence tn. For very small values of 
Up no interesting effect is expected since there are no 
enough pinning centers to stop the vortices. For inter- 
mediate densities, e.g. Up = 0.12, we found the curves 
in Fig. ||-left. The asymptotic value of /^(in) increases 
with Hp. Eventually, if Up is so large that the distance 
between pinning centers is smaller than the range of the 
pinning force, i.e. Up > 7r^^(2ap)^^, there is no differ- 
ence between ZFC and FC preparations and after a very 
few iterations the two curves collapse. 

The maximum intensity in the ramp, /max, is also an 
important parameter. The asymptotic value of /^(in) 
decreases with /max- For a stronger /max, e.g. /max = 1-6 
in Fig. |2|-right, the memory of the initial condition is very 
quickly lost since the vortices are pushed too strongly by 
the force. The perturbation itself, i.e. the sequence of 
current/force pulses, is driving the system to different 
steady states characterized by /^(in — *■ oo). 
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FIG. 3. Vortex trajectories of the FC (left) and ZFC (right) 
samples during the first pulse (lines). Pinning centers where a 
vortex gets pinned are represented with triangles, those from 
which a vortex depins are represented with squares and those 
where vortices get pinned and later depin are represented with 
pentagons. Pinning centers unaffected by the pulse (whether 
there is a vortex permanently pinned to them or not) are 
represented with circles. 

Analyzing the motion of vortices during each pulse 
helps us understanding these results. During the first 
pulse a fraction of vortices in the FC sample remains 
pinned while the rest rapidly depin and move in the 
sample. How many of them depin depends on Tp and 
-^MAx- On the left of Fig. || we display the vortex trajecto- 
ries with lines and the positions of pinning centers with 
symbols; 15 vortices move while the rest remain pinned. 
Some trajectories join two pinning centers. During the 
quiescent time tw a number of vortices may get again 
pinned, depending on Up, the strength of the different in- 
teractions and thermal fluctuations. In the second pulse 
new vortices depin and, consequently, /th decreases. 

Instead, the initially ordered sample rapidly reaches a 
smectic flow (Fig. fright). During this moving pe- 
riod, the vortices that move close to pinning centers get 
pinned; other vortices are pulled by the force and depin, 
but get later pinned by another pinning center. By the 
end of the first ramp, the hexagonal lattice has been de- 
formed and 6 vortices are pinned. At T = the displace- 
ment during the quiescent time is rather small. This 
mechanism repeats in subsequent periods and the thresh- 
old intensity increases in time. This picture is confirmed 
by the study of the position of each vortex relative to its 
closest pinning center as a function of t„. In the initially 
pinned sample, on average and in time, vortices move 
away from the pinning centers. In an initially depinned 
sample, vortices get pinned by a pinning center. A De- 
launay triangulation analysis shows that at long times, 
e.g. tn — 1500, both initially ordered and pinned systems 
become a deformed hexagonal lattice. 

We wish to check if this is a plausible explanation of 
the results in jl^. Certainly, we have almost continu- 
ously stirred the sample while the probing pulses applied 
experimentally, being widely separated in time, are much 
weaker. Unfortunately we cannot match the experimen- 
tal conditions t^ ^ Tp {t^ = 30min = 10^ iterations!). 
On the other hand, mechanical vibrations have not been 
totally eliminated from the experiment and might induce 
an important supplementary means for relaxation fl^ . 
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Another important difference between our numerical 
study and the experiment in [ p^ is that we have worked 
at T = (with the aim of reducing the thermal noise and 
improving the averaging) while experiments are done at 
low but non-vanishing T. The T = protocol is very 
close to the ones used to probe the dynamics of granu- 
lar matter. Indeed, temperature is totally irrelevant in 
these systems and the rearrangement of grains is induced 
by pumping energy in the sample in, e.g., the form of 
periodic taps If the vortex system at T = be- 

haves as granular matter, all the dynamics must then 
take place during the probing pulses while the systems 
should be essentially static in a metastable configuration 
during tw In Fig. ^ we analyse the dependence of the 
relaxation time-scale on by comparing /th'^^'^I^u) for 
t.u, = 3, 144, 288 and Tp — 14.4 in all cases. The longer t^, 
the slower the approach to the asymptotic value. The in- 
set shows /th against the number of perturbing cycles n. 
The three curves now collapse showing that at T = the 
evolution of /th takes place during Tp and it is almost 
completely due to the perturbation. It would be very 
interesting to check if this also hauDens exoerimentallv. 
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FIG. 4. The threshold intensity against time for FC 
(filled symbols) and ZFC (open symbols) preparations and 
tw = 3,144,288, Tp = 14.4 (circles, diamonds, triangles). In 
the inset, the same curves are plotted in function of the num- 
ber of current ramps applied and they rescale perfectly. 

The mechanism for the very slow relaxation in the ex- 
perimental system may be the one put forward in the 
previous paragraph. The relaxation of the FC and ZFC 
samples are controlled by the parameters in the model 
{Up, /max, etc). Hence, by changing the values of the 
parameters, one can explain the difference in the exper- 
imental observations described in the introduction. We 
have also simulated 3D systems with different number 
of layers (1 to 4) and interlayer coupling strengths s (1 
to 9), and observed no qualitative difference in the /^h 
relaxation though its time-scale increases considerably. 
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